Γ(x) [15 D] Alvaro

5123ffa0-03ff-4324-a846-71f4001c3e1d 9 6 12 decimal

&[DATE] &[TIME] - &[FILENAME]

Γ(x) [15 D] Alvaro

Γ.κ 0.99999999999999709182 57.156235665862923517 59.597960355475491248 - 14.136097974741747174 0.49191381609762019978 - 0.33994649984811888699 104 -^* 0.46523628927048575665 104 -^* 0.98374475304879564677 104 -^* - 0.15808870322491248884 103 -^* 0.21026444172410488319 103 -^* - 0.21743961811521264320 103 -^* 0.16431810653676389022 103 -^* - 0.84418223983852743293 104 -^* 0.26190838401581408670 104 -^* - 0.36899182659531622704 105 -^* 151 mat :

ζ Γ z ζ 1 - : zh z 0.5 + : zg zh 607128 / + : zp zg 0.5 zh *^: 14 mat 2 π * sqrt Γ.κ 1 el Γ.κ κ 1 + el z κ + / κ 114 sum + (* zp 2^* e zg -^* eval 21 line :

Valid for Re(ζ)>0

Normalized Hypergeometric ... PDF, CDF

N n p x PDF N 1 + Γ n 1 + Γ * N n + p - 1 + Γ * p 1 + Γ * x 1 + Γ x - n + 1 + Γ * N n + 1 + Γ * x N + p - 1 + Γ * x - p + 1 + Γ * / :

N 60 :

n 40 :

p 10 :

shift 0 :

match the codes [U, p] up 10 :

x pdf N n p x PDF eval :

data U 0 up 0.1 range : i 1 U rows range vy i el N n up U i el PDF eval : for U vy augment 31 line :

Continuous scalar PDF Hypergeometric.

N n 5 x shift - PDF N n 10 x shift - PDF N n 15 x shift - PDF 3 1 sys

data Σ Finite differences Cumul Integral vx data 1 col : vy data 2 col : Y 1 el 0 : i 2 vx rows range Δx i el vx i el vx i 1 - el - : Y i el Y i 1 - el 12 / vy i 1 - el vy i el + (Δx i el * (* + : 21 line for vx Y augment 61 line :

user 51 :

∫user data Σ user row : 5 0.748891 12 mat

hyp N n p ∫user 11 el PDF : 0.2076

Finite differences Cumulative Integral

data Σ data ∫user 1 1 el 1 ∫user 1 1 el 0.25 - 2 2 mat ∫user 1 1 el ∫user 1 2 el + 20 black 1 5 mat ∫user 1 1 el hyp + 20 black 1 5 mat 5 1 sys

Finite Differences [FD] merits . . .

1. ! insensitive ! to noisy data.

2. very accurate in its simplest code.

3. rescues otherwise Simpson reject.

4. recovers 'user' values.

Generally, FD is slightly low. In any

instance, FD needs XY data set only.

]]>

Carlos adaptive integration

]]> ∫Adaptive(f,a,b,N,ε)

f a b iter ε ∫Adaptive sum 0 : i 1 : tol i el ε : aa i el a : h i el b a - (2 / : 51 line .... fa i el a f : fc i el a h i el + f : fb i el b f : s i el h i el fa i el 4 fc i el * + fb i el + (3 / * : L i el 1 : 51 line 13 mat above => initialise the while loop i 0 > fd aa i el h i el 2 / + f : fe aa i el 3 h i el * 2 / + f : 12 mat s1 h i el fa i el 4 fd * + fc i el + 6 / (* : s2 h i el fc i el 4 fe * + fb i el + 6 / (* : 12 mat v1 aa i el : v2 fa i el : v3 fc i el : v4 fb i el : v5 h i el : v6 ε : v7 s i el : v8 L i el : 18 mat i i 1 - : s1 s2 + v7 - abs v6 < sum sum s1 s2 + (+ : v8 iter ≥ Solution Insufficient Number of Intervals : i i 1 + : aa i el v1 v5 + : fa i el v3 : fc i el fe : fb i el v4 : 14 mat h i el v5 2 / : tol i el v6 2 / : s i el s2 : L i el v8 1 + : 14 mat i i 1 + : aa i el v1 : fa i el v2 : fc i el fd : fb i el v3 : h i el h i 1 - el : s i el s1 : 16 mat L i el L i 1 - el : 61 line if if 51 line while sum 41 line eval :

end 5 :

normalized PDF x f 6040 up x PDF :

Program algo style is a bit slower than open style, more convenient. Adapted from Carlos code [20180110]

'N' iterations
precision iter 1000 : ε 1012 -^: 21 sys

User inputs

lower limit [any]
upper limit ≤ end a 0 : b end : 21 sys

t0 1 time :

sanity vs FD f a b iter ε ∫Adaptive 0.748978

1 time t0 - 45.3 s *

Asymptotic sanity ...

end 20 :

normalized PDF x f 6040 end x PDF :

'N' iterations
_ precision _iter 1000 : ε 1012 -^: 21 sys

lower limit [any]
upper limit ≤ end a 0 : b end : 21 sys

t0 1 time :

asymptotic => 1  in cause PDF end flatness f a b iter ε ∫Adaptive 0.999998

1 time t0 - 91.4 s *

an easy test ...

infinite limit x f 1 x / :

'N' iterations
_ precision _iter 1000 : ε 1015 -^: 21 sys

lower limit [any]
_ upper limit _a 104 -^: b 1 : 21 sys

t0 1 time :

sanity 7 D [MCD adaptive] f a b iter ε ∫Adaptive 9.2103404

1 time t0 - 109.5 s *

Simpson [acc 1000] ! 1 x / x 0.0001 1 int 22.7338866

infinite limit x f x Γ :

'N' iterations
_ precision _iter 1000 : ε 1015 -^: 21 sys

lower limit [any]
_ upper limit _a 104 -^: b 1 : 21 sys

t0 1 time :

sanity  8.9637038 [MCD adaptive] f a b iter ε ∫Adaptive 8.9637042

1 time t0 - min 5

infinite limit x f x 0.5 1 -^x - exp * :

'N' iterations
_ precision _iter 1000 : ε 1015 -^: 21 sys

lower limit [any]
_ upper limit _a 1012 -^: b 1 : 21 sys

t0 1 time :

sanity Wolfram ... 1.49364826562485 f a b iter ε ∫Adaptive 1.4936463

1 time t0 - min 2.7

π sqrt 1 erf *

acc dependent examples . . .

]]>

===================================================== Example_1

Romberg integrator companion

f a b Int x f x a b int 32 / b a - (* w G w 01 int * ≡ memo ≡ w s b a - (34 / 2 w * 1 - (* 14 / 2 w * 1 - (3^* - 12 / + (* a + : w G w s f 12 w * 1 - (2^- (* : 32 / b a - (* w G w 01 int * 41 line :

Test function x f x sqrt exp :

7 true digits ... insufficient
integration subspaces. x f x 02 int 5.40744702348009

x f 0 2 1 sys

sanity check the quadruplet Mathcad 11 Romberg 5.40752818465632 ≡ Technical accuracy 5.40752818465627 ≡ Smath Romberg 5.40752818465638 ≡ Smath Integrator 5.4075281 ≡ 41 line

f 02 Int 5.40752818690417

===================================================== Example_2

"Romberg singular" ... option integral accuracy: 10000 f a b Int w s b a - (34 / 2 w * 1 - (* 14 / 2 w * 1 - (3^* - 12 / + (* a + : w G w s f 12 w * 1 - (2^- (* : 32 / b a - (* w G w ε 1 int * 31 line :

ε 1015 -^:

x f x 1 - x ln + (* :

t0 1 time :

Smath Romberg singular end point f a e Int 1.84726399729247 -

1 time t0 - min 0

Mathcad Adaptive 1.84726402473266 - ≡

===================================================== Example_3

x f 41 x 2^+ / :

x f

NOT a case for Romberg. The integrand
is perfectly suited for Simpson. π x f x 01 int f 01 Int 31 mat 3.14159265358979 3.14159265358979 3.14159265283962 31 mat

===================================================== Example_4

x f 2 x * cos x exp / :

x c 10 x ≤ (x 2 π * ≤ (& cases :

x f x c * 2 π * 0.5 2 π * 0.5 - 2 2 mat 2 1 sys

x f x 02 π * int f 02 π * Int 21 mat 0.199626452011026 0.199626526204812 21 mat

===================================================== Example_5

x f x - exp x / :

x c 11 x ≤ (x 100 ≤ (& cases :

x f x c *

Integrand suitable for Romberg transform. x f x 1100 int f 1100 Int 0.219383934395520 31 mat 0.22071060184138 0.219386700436576 0.21938393439552 31 mat

===================================================== Example_6

x f 11 x 2^+ / :

x c 10 x ≤ (x 3 ≤ (& cases :

x f x c *

x f x 03 int f 03 Int 21 mat 1.24904577238205 1.249045772194 21 mat

three exact Simpson . . .

]]>

u b 0 u 2 - ≤ 16 / 2 u + (3^* 2 - u ≤ (u 1 - ≤ (& 16 / 2 u + (3^* 23 / 1 u + (3^* - 1 - u ≤ (u 0 ≤ (& 23 / 1 u - (3^* - 16 / 2 u - (3^* + 0 u ≤ (u 1 ≤ (& 16 / 2 u - (3^* 1 u ≤ (u 2 ≤ (& 0 cases :

The Gaussian b-spline

Gaussian basis 'b-spline'

2 - x ≤ x 2 ≤ & x b Undefined if

exact acc 100 x b x 2 - 2 int 1

rk cumulative integrator . . .

]]>

x y D x pdf :

sol 0010100 x y D rkfixed :

x pdf sol 2 1 sys

arc-length Engineering . . .

This example applies to 2D b_spline, a most

common polyline. The data vector can be obtained

from more common 1D spline wrt the project.

Offered FD [Finite Differences], thus universal

for any type of the function to be arc-measured.

]]> plotG

vx vy a b c plotG a char ≡ b size ≡ c color ≡ 13 mat n vx length : plot vx 1 el vy 1 el a b c augment : i 2 n range plot plot vx i el vy i el a b c augment stack : for plot 51 line :

u b 0 u 2 - ≤ 16 / 2 u + (3^* 2 - u ≤ (u 1 - ≤ (& 16 / 2 u + (3^* 23 / 1 u + (3^* - 1 - u ≤ (u 0 ≤ (& 23 / 1 u - (3^* - 16 / 2 u - (3^* + 0 u ≤ (u 1 ≤ (& 16 / 2 u - (3^* 1 u ≤ (u 2 ≤ (& 0 cases :

The Gaussian b-spline

Gaussian basis 'b-spline'

2 - x ≤ x 2 ≤ & x b Undefined if

0 b 0.6667

u v Bsp v m el u m - b * (m 1 v rows sum :

observe u b u 2 - 2 int 1

X 1.875 1.875 1.95 2 2.05 2.125 2.125 2.125 18 mat : Y 11101 - 1 - 1 - 1 - 18 mat : 21 line

The XY plot allows to plot in the reverse X direction. The choice of the control points is by "user". You may have to trick the ranging of 't'. "data" is an export. There is no limit in the length of X, Y. With this kind of spline, you always have to at least duplicate end points [here RHS is triplicate]

XY X transpose Y transpose augment : u oy u Y transpose Bsp : u ox u X transpose Bsp : 31 sys

Points X transpose Y transpose . 14 black plotG :

data t 1.75 6.75 1.8 range : i 1 t rows range vx i el t i el ox : vy i el t i el oy : 21 line for vx vy augment 31 line :

data rows 101

XY data 2 1 sys XY data Points 3 1 sys

Finite differences derivative X data 1 col : Y data 2 col : 12 mat n X rows : deriv 1 el 0 : deriv n el 0 : i 2 n 1 - (range deriv i el 0.5 Y i 1 - el Y i el - X i 1 - el X i el - / Y i el Y i 1 - el - X i el X i 1 - el - / + (* : for deriv 1 el deriv 2 el : deriv n el deriv deriv rows 1 - el : D X deriv augment : Arc integrand .................... xx D 1 col : δx D 2 col : 12 mat i 1 δx rows range yy i el 1 δx i el (2^+ sqrt eval : 11 line for Arc xx yy augment : Arc length ....................... vx Arc 1 col : vy Arc 2 col : 12 mat S 1 el 0 : i 2 vx rows range Δx i el vx i el vx i 1 - el - : S i el S i 1 - el 12 / vy i 1 - el vy i el + (Δx i el * (* + : 21 line for Σ vx S augment : 181 line

Derivative of the function

Arc integrand

Arc

Cumulative arc length

Because

of the ranging of 't' that is difficult

to adjust, spooling to file is a

desirable conservative option.

]]>

CWD C:\SmathFile CurrentDirectory :

Arc CWD ArcRobotics wfile 1

Arc CWD ArcRobotics rfile :

X Σ 1 col : Y Σ 2 col : 12 mat

arc length from interpolation x arc X Y x ainterp :

2 arc 1.0055

x f 41 x 2^+ / :

x f

NOT a case for Romberg. The integrand
is perfectly suited for Simpson. π x f x 01 int f 01 Int 31 mat 3.14159265358979 3.14159265358979 3.14159265283962 31 mat